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ABSTRACT 

In determining Mercury's core structure from its rotational properties, the 
value of the normalized moment of inertia, C/MR"^, from the location of Cassini 1 
is crucial. If Mercury's spin axis occupies Cassini state 1, its position defines the 
location of the state, where the axis is fixed in the frame precessing with the orbit. 
Although tidal and core-mantle dissipation drive the spin to the Cassini state with 
a time scale O(IO^) years, the spin might still be displaced from the Cassini state 
if the variations in the orbital elements induced by planetary perturbations, which 
change the position of the Cassini state, cause the spin to lag behind as it attempts 
to follow the state. After being brought to the state by dissipative processes, the 
spin axis is expected to follow the Cassini state for orbit variations with time 
scales long compared to the 1000 year precession period of the spin about the 
Cassini state because the solid angle swept out by the spin axis as it precesses 
is an adiabatic invariant. Short period variations in the orbital elements of small 
amplitude should cause displacements that are commensurate with the amplitudes 
of the short period terms. The exception would be if there are forcing terms in 
the perturbations that are nearly resonant with the 1000 year precession period. 
The precision of the radar and eventual spacecraft measurements of the position 
of Mercury's spin axis warrants a check on the likely proximity of the spin axis to 
the Cassini state. How confident should we be that the spin axis position defines 
the Cassini state sufficiently well for a precise determination of C /MR^l 

By following simultaneously the spin position and the Cassini state position 
during long time scale orbital variations over past 3 million years (Quinn et ai, 
1991) and short time scale variations for 20000 years (JPL Ephemeris DE 408, 
E. M. Standish, private communication, 2005), we show that the spin axis will 
remain within one arcsec of the Cassini state after it is brought there by dissipa- 
tive torques. In this process the spin is located in the orbit frame of reference, 
which in turn is referenced to the inertial ecliptic plane of J2000. There are no 
perturbations with periods resonant with the precession period that could cause 
large separations. We thus expect Mercury's spin to occupy Cassini state 1 well 
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within the uncertainties for both radar and spacecraft measurements, with corre- 
spondingly tight constraints on C /ME? and the extent of Mercury's molten core. 
Two unlikely caveats for this conclusion are 1. an excitation of a free spin preces- 
sion by an unknown mechanism or 2. a displacement by a dissipative core mantle 
interaction that exceeds the measurement uncertainties. 



1. Introduction 

Radar observations have begun the process of determining the extent of Mercury's molten 
core (Margot et al. 2003), and the experiment will be completed when the MESSENGER 
spacecraft orbits Mercury in 2011 (Solomon et al. 2001). The BepiColombo spacecraft will 
complement and augment the observations of MESSENGER, but will orbit Mercury sometime 
after 2011 (Anselmi and Scoon, 2001). The experiment is based on the product of three factors 
(Peale, 1976; 1981; 1988; 2005; Peale et al. 2002). 

\B-A)\ C )\ MR^ J C - ' ^ ' 

where A < B < C are the principal moments of inertia of Mercury with being the 
polar moment of inertia of the mantle alone, and M and R are Mercury's mass and radius 
respectively. The first factor is determined by the amplitude of the physical libration, 0o — 
[3{B — A)/2Cm]{l — lle^ -I- • • •), where e is the orbital eccentricity and Cm appears in the 
denominator because the liquid core is not expected to follow the mantle during the short 
period librations (Peale, et al. 2002). The second factor follows from the analysis of generalized 
Cassini's laws for Mercury (Colombo, 1966; Peale, 1969; Beletskii, 1972) (see Section 3). 

C ^ nJ2 /(e) 
Mi?2 (sin/')/ic + cos/'' ^ ^ 

where /(e) = (^210(6) + 2C22G2oi(e) / J2, n — -^/GM^Jc^ is the orbital mean motion, (G is the 
gravitational constant, Mq is the solar mass, a is the semimajor axis of the orbit), the functions 
G'2io(e) = (1 - e2)-'V2 and G'2oi(e) = 7e/2 - 123eVl6 + ... are defined by Kaula (1966), J2 
and C22 are the second degree coefficients in the expansion of Mercury's gravitational field, /' 
is the inclination of the orbit plane to the Laplace plane, the plane on which Mercury's orbit 
precesses with nearly constant inclination and at a nearly constant rate whose magnitude is 
and ic is the obliquity of the Cassini state. The orientation of the Laplace plane and the 
value of are determined by averaging the orbital element variations over a suitable time 
interval as discussed below. The last factor is found from C22 — {B — A)/{AMR^). 

Fig. 1 shows the geometry of Cassini state 1 for Mercury, where Mercury's obliquity i = ic 
with the ascending nodes of the equator plane on the orbit plane and of the orbit plane on the 
Laplace plane remaining coincident as the spin and orbit normal precess around the normal to 
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the Laplace plane. The Laplace plane is determined mainly by Venus, Earth and Jupiter as is 

the precession of Mercury's orbit. Because of the orbit precession, the spin precesses around 
the Cassini state and tends toward that state from dissipative effects rather than toward the 
orbit normal (Peale, 1974: Ward, 1975). 




Fig. 1. — Geometry of Cassini state 1. The ascending node of the equator on the orbit 
plane and the ascending node of the orbit on the Laplace plane remain coincident as they 
precess around the normal to the Laplace plane at the current rate corresponding to the orbit 
precession period of about 300,000 years. 



The determination of 00, <^2, (^22, and ic thus yields Cr„/C. Empirical constraints on C/Mi?^ 
must wait for MESSENGER determinations while in orbit. Grucial to the determination of 
moment of inertia is the assumption that the spin axis occupies Gassini state 1 (Eq. (2)). 
Only if we are confident that the measured obliquity i — ic, can we constrain the value of 
C/MR^ with sufficient precision to determine the extent of Mercury's liquid core. Knowing the 
internal structure of Mercury is important for inferring the thermal history, overall chemical 
composition and the constraints Mercury can place on the details of origin and evolution of 
the terrestrial planets (Harder and Schubert, 2001; Solomon et al. 2001). 
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Tidal friction will drive Mercury's spin to Cassini state 1 from virtually any initial obliquity 
(Peale, 1974; Ward, 1975), and the addition of dissipation between a liquid core and solid 
mantle will hasten that evolution, provided the initial obliquity i < 90°, with a time scale 
from both tides and core-mantle dissipation together of order 10^ years (Peale, 2005). As 
Mercury's spin axis approaches the Cassini state, it will precess around that state with a 
period near 500 or 1000 years depending on whether or not the liquid core is dragged along 
with the spin precession. A finite amplitude free precession would frustrate the determination 
oi ic, but the short time scale for its decay makes a remnant free spin precession unlikely, 
unless there is a recent, unspecified excitation mechanism (Peale, 2005). 

Another possibility for a displacement of the spin from the Cassini state is if the spin is 
unable to follow the Cassini state sufficiently closely as the position of the state changes due to 
the variations in the orbital elements and Laplace plane orientation from planetary perturba- 
tions and changing solar system geometry. We can expect the spin to remain reasonably close 
to the Cassini state after being brought there by dissipative processes because the following 
action integral is an adiabatic invariant, where the integral is nearly constant if the period of 
the angle variable is short compared to all other relevant time scales (Goldreich and Toomre, 



where the action variable p' = —a{l — cos 6') with a being the spin angular momentum, and 
q' — (f)' is the angle variable, with 9' and 0' being the ordinary spherical polar coordinates of 
the spin vector in the system with the Z axis aligned with the Cassini state. That p' and q' are 
conjugate variables is verified by dp' /dt — —dH/dq' and dq'/dt — dH/dp' , with H being the 
Hamiltonian for the Mercury's rotational dynamics written in the frame precessing with the 
orbit (Peale, 1974). The action integral is seen to be —a times the solid angle contained by the 
spin vector as it precesses around the Cassini state, and it is approximately conserved if the 
spin precession rate 0' is fast relative to the significant changes in the parameters determining 
the position of the state. The angle variable describes the precession of the spin about the 
Cassini state. We shall see below that the large amplitude variations in the orbital parameters 
relative to the ecliptic plane have periods exceeding 5 x 10^ years, which is sufficiently long 
compared to the 1000 year precession period that one expects the adiabatic invariant to keep 
the spin close to the current position of the Cassini state (once it is there) as the latter's 
position changes slowly on these time scales. The adiabatic invariant is not conserved on 
the time scale of short period variations, but these variations are of small amplitude, and we 
shall sec that the spin follows the Cassini state defined by the orbital elements averaged in a 
2000 year window over the 20,000 year JPL Ephemeris DE 408. The precision of the radar 
determinations of Mercury's spin properties, and that anticipated for the MESSENGER and 
BepiColombo missions warrants a check on just how well the spin axis follows the Cassini 
state for all variations of the parameters defining the state. 



1969; Peale, 1974). 




(3) 
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Our purpose here is to develop a formalism to test how close to the Cassini state the spin 
axis remains for both slow and fast variations in the orbital parameters. We shall consider 
the variations of e, /, Q, and their derivatives due to the planetary perturbations, where / 
and Q are the are the orbit inclination and longitude of the ascending node of the orbit on 
the ecliptic plane respectively. Dissipative processes will relentlessly drive the spin to the 
Cassini state, so the test will determine how closely Mercury's spin axis follows the Cassini 
state position as the position of the state changes due to the orbital element variations. To 
this end we develop the equations of motion of Mercury's spin vector relative to the orbit 
frame of reference, averaged over an orbit period, in Section 2, and find the position of the 
spin and the position of the Cassini state as a function of time as variations in the orbital 
parameters ensue. Finding the Cassini state position will involve defining a set of coplanar 
vectors that includes the normal to the Laplace plane. The inclination of the orbit relative to 
the Laplace plane and the rate of precession of the ascending node of the orbit on the Laplace 
plane do not have to be known with high precision to define the position of the Cassini state 
accurately (Yseboodt and Margot, 2005). For the slow variations in e, /, fl, dl/dt and dfl/dt, 
we shall vary the orbital elements and their time derivatives according to simulations by T. 
Quinn (Quinn et al. 1991) over the past 3 x 10^ years. ^. For the short period variations in the 
orbital elements we follow the spin and Cassini state positions as e, /, Q, dl/dt and dD,/dt 
vary according to the 20,000 year JPL Ephemeris DE 408 provided by Myles Standish, but 
now the Cassini state position will be determined by elements and rates averaged over the 
2000 year window mentioned above. Periodic variations of the same variables with amplitudes 
and periods representative of the true variations are also applied to identify the sources of the 
fiuctuations in the spin-Cassini state separation that are observed. 

In Section 3 we determine the time varying position of the Cassini state, where variations in 
dl/dt and dfl/dt change the Laplace plane orientation and the variations in /, e and w change 
ic- This determination allows us to compare the position of the spin vector with the Cassini 
state position at any time. All dissipative forces are ignored, and principal axis rotation is 
assumed. The neglect of dissipation means we couple the liquid core firmly to the mantle, 
which is equivalent to Mercury's having a solid core. Consequences of the relaxation of this 
latter assumption will be pursued in a later paper. 

In Section 4, the spin vector, placed initially in the Cassini state, is shown to remain within 
approximately 1" of the Cassini state position as this position is continuously redefined for 
the continuously varying Laplace plane and orbital elements as given over the past 3 x 10^ 
years by the Quinn simulation. Equivalently, the spin vector would maintain any small initial 
separation from the Cassini state to within 1". The effect of higher frequency terms on the 
spin-Cassini state separation is found in Section 5. Large separations of the spin from the 
Cassini state are forced for periodic variations in the inclination with periods near the spin 



^ftp://ftp. astro.washington.edu/pub/hpcc/QTD 
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precession period. However, any terms with periods near resonance with the spin precession 
period are not evident in the real variations of e, /, fl, dl/dt and dfl/dt . This is demonstrated 
by the spin and Cassini state positions (from locally averaged values of elements and rates) 
remaining within 1" for variations of these parameters given by JPL Ephemeris DE 408. A 
summary and discussion follows in Section 6, where we end with a conjecture on the possible 
effect of the liquid core on the position of the spin axis. 



2. Equations of variation 

Fig. 2 defines the coordinate systems centered on Mercury to be used along with some 
of the variables. The XYZ system is the inertial echptic plane of J2000 with the Z axis 
perpendicular to the plane. The X'Y'Z' system has the Z' axis perpendicular to the orbit 
plane and the X' axis along the ascending node of the orbit plane on the XY plane. The xyz 
system is the principal axis system fixed in the body. The orbit plane is inclined to the XY 
plane by angle /, Q is the longitude of the ascending node of the orbit plane on the XY plane 
measured from the inertial X axis along the direction to the vernal equinox, Gq = e^' is the 
orbit normal, = is a unit vector along the spin axis, fls is the longitude of the ascending 
node of Mercury's equator on the orbit plane relative to the X' axis, i is inchnation of the 
equator plane relative to the orbit plane (obliquity), ip is the angle between the ascending 
node of the equator and the x principal axis, r points toward the Sun in the orbit plane and 
cu and / are the argument of perihelion and true anomaly locating the Sun in the orbit plane 
relative to the X' axis. 

In Peale (2005) it is shown that although the actual spin precession trajectory is slightly 
elliptical and the rate of precession is not quite uniform, a good approximation for the time 
variation of the spin vector is given by 

des 

—T^ — Ki cos i {bs X eo) = —Ki cos i sin i(cos fts e^' + sin Jl^ e^') (4) 

CLv 



where 



^J2G2io(e) + 3C22G2oi(e) 



Kf{e), (5) 



with K = nMR^J^/C and /(e) = ^210(6) + 2(722^201 (e)/J2, (used in Eq. (2)) and with 
■0 = 3n/2 being assumed. Eq. (4) follows from the average of the equations of motion over 
the orbit period as mentioned above with ip — 1.5n, and neglect of terms with coefficients 
that are small compared with Eq. (5). We can consider only the variation of the unit vector 
because ip = 1.5n from the spin orbit resonance if we neglect the small physical librations. 
Eq. (4) shows the expected regression of the spin vector about the orbit normal. It is easiest 
to determine the equations of motion of the spin in the orbit frame of reference starting from 
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Fig. 2. — Coordinate systems used in determining time variation of Mercury's spin vector, 
where and eg are unit vectors along the spin axis and orbit normal respectively. The XY 
inertial plane is the echptic plane of J2000. 

the relation ^ ^ 

— ^ = — ^ + w X Bs ^ Kicosi {eg X Bo), (6) 

dt inertial Ctt orbit 

where = sin i sin fl^ ^x' — sin i cos fls e^' + cos i Bz' and Bo = Bz' are defined in Fig. 2. The 
Bi are unit vectors along the respective axes. 

The angular velocity w is nominally w/j,, the angular velocity of the orbit normal about the 
Laplace plane normal, but we can make another choice that simplifies the equations. The 
precession of the orbit normal about the Laplace plane normal leads to a velocity of the 
orbit normal Bg of v — W]^ x b^, where is the angular velocity of the orbit normal parallel 
to the Laplace plane normal. Differentiation of the components of Bg in the ecliptic frame 
yields components of v that can be set equal to the like components of the cross product. 
However, these equations for the components of are not linearly independent, so one can 
only solve for two of the components in terms of the other, which serves as a free parameter. 
A convenient choice of the free parameter is the Z component of w, where we have abandoned 
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the subscript, since another constraint is necessary to determine the Laplace plane normal 



w 



dl ^ / dQ. r ■ 

— COS i I + [Wz — -7-) tan 1 sm\l 

CLL CLL 



— smO + {—Wz + -7-) tan/ cos y 

CLL CLL 



BY + Wzez- 

(7) 



All of the solutions for w lie in the plane determined by the orbit normal and the Laplace 
plane normal, where the latter is one of the set. The coplanarity is understood since all the 
w X eg must produce the same v. If /' is the angle between and w, v — w sin/' is a 
constant, and the magnitude w must decrease as /' increases. Eq. (7) is Eq. (13) of Yseboodt 
and Margot (2005), where they determine tu^ at a given epoch by adding the numerically 
determined constraint that the variation in the inclination /' is minimized in a 2000 year 
window centered on the epoch for data obtained from the 20,000 year JPL Ephemeris DE 408. 
At the epoch J2000, the ecliptic latitude and longitude of are 86.725° and 66.6° respectively 
corresponding to Wz = —1.91 x 10^^ radians/year, /' = 8.6° and an instantaneous precession 
period of 328,000 years. There is some uncertainty in these values from the statistical nature 
of the minimization process, but we shall see below that the uncertainty in the Laplace plane 
normal does not compromise the determination of the Cassini state. 

The simplest form for w in Eq. (7) is for Wz — dfl/dt. Since any choice of Wz in Eq. (7) 
yields the correct instantaneous motion of Cq, we can use this simplest form, 

^ dl . ,dVL dVL 

w = -776^' + sm /— ey' + cos I-rr^z', (8) 
dt dt dt 

in Eq. (6) to determine the equations of motion of the spin in the orbit frame of reference. 
This choice of expressed here in the orbit frame, is just the vector sum of dl /dt and dD,/dt 
(Fig. 3). Equating like components in Eq. (6) yields three linearly dependent equations that 
can be solved uniquely for di/dt and dVL^/dt. There results 

di ■ T ■ dVl dl 

— — — smismiZB— cosUe — 

dt dt dt 

dfls cos i sin dl 

—— = -/ricosz + — — — — 

dt sm I dt 

cos i cos sin / + sin i cos / dO. 

sin i dt ' 

where dl /dt and dil/dt are assumed known. 

As expected, di/dt = and dfls/dt — —Kicosi — constant if dl/dt and dVl^^/dt are both 
zero. If dl/dt — but dVt/dt is a negative constant, the equations describe the motion of 
the spin precession about a Cassini state in a uniformly precessing orbit. The Cassini state is 
displaced from the orbit normal in the plane defined by the orbit normal Gq and the normal 
= —Wi^/wi^ (Fig- 1) (Colombo 1966, Peale, 1969). This is the same plane containing all of 
the vectors defined by Eq. (7). 
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Fig. 3. — The vector sum of (Kt/dt and dl/dt yields an effective precessional angular velocity 
w of the orbit normal, which is the simplest choice from Eq. (7). 

To eliminate the sin i singularity in Eqs. (9) we introduce the variables p — sini sin Jl^ and 
q — sini cos Jl^, with the result 

— — i^/(e)g-\/r^-p^^-^ — (sin /-y/r^-p^^-^ + gcos /)^, 
^ = Kf{e)p\J\ - 'p^ - (f - \J\- - (f'^ +pcos/^. (10) 

The numerical solution of Eqs. (10) yields the position of the spin vector of Mercury in the 
X'Y'Z' orbit system when all of e, / and Q are varying. 



3. Cassini state position 



We wish to start the spin axis either in the initial Cassini state or close to it and determine 
how close it remains to the changing position of the state as the spin varies according to Eqs. 
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(10). Prom Eq. (12) of Peale (1974), the position of the Cassini state 1 is defined by 



^l-el^, cos /' + e,,, sin /' - 2i?'ecz' ^l-e^^, - 25(1 + e^^,) ^l-e^^, 





0, (11) 



where 



^ - -I—P, ] -'2Lr210(ej, 

^ - l^^C,,GMe), (12) 

with Bcx'.y'.z' being the components in the orbit system {X'Y'Z') of Be, the unit vector along 
the Cassini state. Por Mercury, Cassini state 1 is very close to the orbit normal {e.g., Peale, 
1969, 1974), which geometry is shown in Pig. 1. We can then set ^/l — e^^, = | sin id ~ ic, 
Ccz' ~ 1 to first order in ic and write 

w sin /' 



w(2i?' + 45 - cos /')' 



(13) 



where we have replaced w^hy w — •>/ {dl/dty + {dQ/dty, the magnitude of the vector defined 
in Eq. (8) and used in the equations of motion of the spin. The inchnation F in Eq. (13) 
defined below is distinct for each choice of w. But since we have shown above that wsinl' 
has the same value for all the values of w compatible with Eq. (7), the only change remaining 
after this substitution is the change in the term w cos /' in the denominator. Wc shall sec 
below that 2R' + 4S > 250 ^ cos/', so there is less than ~ 0.1% change in (~ 0.05" 
at the current epoch) effected by this substitution. Like the equations of motion for the spin 
(Eqs. (10)), use of this w in place of greatly simplifies the numerical determination of the 
Cassini state position. 



The inclination of the orbit /' to use in Eq. (13) is defined by 

— (co • ty) cos/ 

w ~ ^/i + {di/dt)y{dn/dty 



cos/', (14) 



where the minus sign on the left hand side defines /' as a small positive angle when dD,/dt < 0. 
Eq. (13) thus defines the obhquity of the Cassini state in the orbit frame of reference for the 
values of e, /, dl/dt and dD,/dt at a particular time. As noted above, the Cassini state is in 
the plane defined by Bq = e^' and w, which plane also contains Wi^. 

Using Eq. (8), we can write (Fig. 2) 

dt dt ^ ^ 

to define the components of a vector in the orbit [X'Y'Z') system that is perpendicular to the 
plane containing the instantaneous Cassini state 1 and w appropriate to the instantaneous 



position of the latter vector. We already have the obhquity of the Cassini state from Eq. (13), 
so a unit vector along the instantaneous Cassini state in the X'Y'Z' system is 



Be = sin ic sin Vt^ e^' — sin cos Vt^ By/ + cos e^' , 



(16) 



where 



cos fic 



-{dVt/dt) sin/ 



sin Q,c 



dl/dt 



(17) 



yJ{dQ/dtYBin'l+{dI/dtY' 



with being the angle between the X' axis and the vector in the X'Y' plane defined by Eq. 



We are now in a position to check just how well the adiabatic invariant is satisfied, or more 
directly, how well the spin follows the Cassini state. In Fig. 4 the variations of Mercury's 
eccentricity and inclination to the ecliptic of J2000 are shown for the last 3 x 10^ years from 
the data kindly provided by T. Quinn (Quinn, et al. 1991, See footnote ^.). These data 
have been filtered to exclude any terms with periods less than 2000 years. The justification 
for this filtering is that such terms usually are small amplitude oscillations that average to 
yield negligible contributions to the element variations. Note that keeping only contributions 
from terms with periods longer than 2000 years excludes any terms close to the spin precession 
period of between 1000 and 1100 years. The consequences of near resonant terms arc discussed 
in Section 5. For the integrations, we divide Eqs. (10) by A', defined after Eq. (5), such that 
Kt is a dimensionless time. The period 2tx / K ^ 1365 years would be close to the spin 
precession period for small e if C22 = 0. The dimensionless equations that are integrated 
are then Eqs. (10) with K removed from the coefficient of the first term on the rhs of the 
equations, and Kt — > t. 

We note that the spin vector 



and that Eq. (16) can be written similarly with p and q replaced by Pc = sin?',. sin and 
Qc = sin cos fie with and being defined in Eqs. (13) and (17) respectively. We wish to 
determine the angle 6 between and Cc as a function of time, where 



(15). 



4. 



Results for slow variations 



65 = sini sin e^' — sin^ cos Gy' + cosi e^'. 




(18) 



cos 5 — eg • ec 




(19) 
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t (yearsxlO"^) 



Fig. 4. — Variation of Mercury's eccentricity and orbital inclination to the ecliptic of J2000 
from T. Quinn (ftp://ftp.astro.washington.edu/pub/hpcc/QTD). 

so that 

5 ^ + ^2 + p2 + q2 _ 2pp^ - 2qq, (20) 

gives a more accurate value for 6 than numerically calculating cos^^ 6 when 6 is very small. 
A cubic spline (Press et al. 1986) through the Quinn data points is used to determine e(t), 
I{t), Q, dl/dt and dfl/dt in the calculation. The fact that contributions to the Quinn data 
with periods less than 2000 years have been eliminated, means the values of the variables can 
be used directly in calculating the position of the Cassini state without additional averaging 
(See Section 5). With the exception of the starting and ending values used to evaluate the 
second derivatives in the spline fit, the values of dl/dt and dfl/dt are determined as defined 
in the spline algorithm. 

We show in Fig. 5 the trajectories of the projections of the unit spin vector and the unit 
vector in the direction of the Cassini state ec on the orbit plane over the 3 x 10^ year interval 
covered by the Quinn data. Initially, the two unit vectors are coincident {6{t = 0) = 6o = 
cos~^(es-ec) = 0). That the two trajectories are almost indistinguishable throughout the time 
interval shows how well the spin axis follows the Cassini state for variations in the parameters 
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Fig. 5. — Trajectories of the projections of the unit spin vector and the unit vector in the 
Cassini state direction Bc on the orbit X'Y' plane for variations in e, /, Q, dl/dt and dQ/dt 
according to 3 x 10^ year simulations by T. Quinn. 

that define the state. More precisely, we show in Fig. 6 the actual variations in 6{t) for 
initial angular separations of 0" and 10". The results for panel a were obtained with spline fits 
through all 6185 data points from the Quinn simulation, from which e{t), I{t), Q{t) dl{t)/dt 
and dQ{t)/dt were determined at each call to the Bulirsch-Stoer integrator and at points 
therein. For panel b a spline fit through only every other data point was used including the 
first and last point of the full data set. For the same initial conditions, the deviations of the 
spin axis from the Cassini state in panel a are more scattered than those in panel b and reach 
a maximum of 2.8", whereas the maximum deviation in panel b is only 1.2", which difference 
is explained as follows. 

The interval between data points in the Quinn simulation is 180,000 Julian days or approxi- 
mately 493 years. The shortest period covered in a Fast Fourier Transform (FFT) of the data 
is the Nyquist period of 986 years (e.g. Press et al. 1986), which is less than the 1000-1100 
year spin precession period. Although the spectral power at periods less than 2000 years are 
suppressed by as much as 9 orders of magnitude in the Quinn data (Quinn et al. 1991), we 
show in Section 5 that the system is extremely sensitive to variations in the inclination at 
periods near the spin precession period, which period is included in the full data set. An FFT 
of dl/dt constructed from the spline fit to the full data set sampled every 200 years, with 
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either Hanning or Welch window functions, yields some power down to the 984 year Nyquist 
period, which includes the precession period. 

If we select only every other data point from the simulation, the Nyquist period is increased 
to 1972 years, which now excludes the resonant period. The FFTs of the spline fit and its 
derivative through the inclination points of the half data set sampled at 200 years reproduces 
the spectral power distribution of the full data set, but drops by more than two orders of 
magnitude for periods shorter than the 1972 year Nyquist period of the half data set. In 
other words, the half data set contains more than 100 times less power at the resonant period 
than the full data set, and this accounts for the more erratic behavior and larger maximum 
separation in panel a. We show in Section 5 that there is little or no spectral power near a 
period of 1000 to 1100 years in the full Mercury cphcmcris. So the separations of the spin from 
the Cassini state indicated in panel 6, where the spin is initially in the Cassini state, is more 
representative of the true limits on this separation as the orbital elements vary according to 
the Quinn solar system simulation over the 3 million year interval. Panel c shows that the spin 
remains within about 1" of an initial separation of 10" indicating that the adiabatic invariant 
is reasonably well preserved. 

The source of the maximum separation of the spin from the Cassini state in panel h of Fig. 6 
can be inferred by holding dVl/dt constant with precession period of 300,000 years and varying 
/ and e according to 

I = 7o + ^sin(^^^, 

e = eo + Aesin (^^ + 0^, (21) 

and use the Quinn data to assign mean values Iq = 6.0°, Cq = 0.19 and the following amplitudes 
with associated periods: 



A, 


= 5.5° 


for 


Pi 


= 1 X 10^ 


yr, 


A, 


= 1.5° 


for 


Pi 


= 2 X 10^ 




A, 


= 0.6° 


for 


Pi 


= 5 X 10^ 


yr, 


Ae 


= 0.06 


for 


Pe 


= 8 X 10^ 


yr, 


Ae 


= 0.01 


for 


Pe 


= 5 X 10^ 


yr- 



(22) 



Fig. 7 shows the results for 6{t) for e and / variations with several combinations of the 
above amplitudes and periods. The initial value (5o = is chosen for all cases. The large 
amplitude, long period variations lead to maximum deviations of the spin from the Cassini 
state of about 0.25", whereas shorter period variations (P = 5 x 10^ yr) with smaller amplitude 
lead to a maximum deviation a little more than 1". It is these shorter period variations in the 
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Fig. 6. — Separation S of the spin vector from Cassini state 1 for variations in e, /, dl/dt and 
dVt/dt from a simulation by T. Quinn. The spin is initially at the Cassini state in panels a 
and b and separated initially by 10" in panel c. The separation 5{t) in a is constructed from a 
spline fit through all 6185 data points in the Quinn simulation, whereas that in h results from 
a sphne fit through every other point or 3093 points total. The reason for the differences in 
a and h is discussed in the text. Panel c shows the constancy of 5 to within ±1" for an initial 
value of 10". 

the orbital parameters that are causing the maximum deviations of the spin from the Cassini 
state in Fig. 6. The fact that different values of with Pg = 5 x 10^ yr in combination with 
Ai = 1.5° with Pf = 2 X 10^ yr yield comparable maximum values of 6{t) shows that variations 
in I are much more important than the variations in e in causing to not follow the Cassini 
state. Next, the maximum deviation from the Cassini state is almost proportional to Aj for a 
fixed Pi, from which we infer that it is really the maximum dl/dt that determines how much 
Bs deviates from Sc- If we set I — Iq + {dl /dt)maxt, with {dl/dt)max — 1-2 x 10~^°/yr from 
the spline fit to the Quinn data, the excursions of the spin away from the Cassini state over 
approximately 3.5 x 10^ years are shown in Fig. 8. For e = 0.19, the maximum 5 < 0.6". 
Imposing a periodic variation in e at the maximum amplitude of 0.06 while / is increasing 
linearly, induces a modulation in the maximum S, whose peaks approach only to 0.8". 
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Fig. 7. — Separation of spin and Cassini state for periodic variations I = Iq + Aj sin (27rt/Pf), 
e = eo + Ae sin {27it/Pe + (f)), where initially the spin axis and Cassini state are coincident. The 
amplitudes and periods are characteristic of those in the Quinn data. Iq = 6.0° and Cq = 0.19 
are the real approximate mean values. These mean values and (p — are assumed for all 
the cases, a) A^ = 5.5°, = 10^ jr, A^ = 0.06, Pe = 8 x 10^ yr; b) Aj = 1.5°, Pj = 2.0 x 
lO^yr, Ae = 0.06, Pe = 8 x lO^yr; c) A^ = 1.5°, = 2 x lO^yr, A^ = 0.01, Pe = 5 x lO^r; 
d) Aj = 0.6°, P, = 5 X lO^r, A = 0.01, Pe = 5 x lO^r- 



5. Results for short period veiriations 



The Quinn data has been filtered to eliminate all variations with periods less than 2000 
years, so it is instructive to look at the response of the system to short period fiuctuations. 
For this purpose we use the complete Mercury orbital element variations, sampled every 500 
Julian days from the JPL Ephemeris DE 408, that were kindly provided by Myles Standish. 
First we determine the location of the current Cassini state. Fig. 9 shows the variation 
over 20,000 years of e, / and Q relative to the echptic plane, of J2000, where is measured 
from the vernal equinox. The variations are dominated by nearly linear secular changes 
with short period fluctuations superposed whose amplitudes are within the line widths of the 
curves. The current position of the Cassini state is necessary for the interpretation of the 
radar and future spacecraft information. We can determine a preliminary position of the 
Cassini state that will be refined when the MESSENGER spacecraft orbits Mercury. The 
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Fig. 8. — Separation of spin and Cassini state for the maximum dl/dt in the sphne fit to the 
Quinn data. / = 7.35° + (1.2 x 10-^° /yr)t. a) e = 0.19. b) e = 0.19 + 0.06 sin 27rt/10V- c) 
e = 0.19 + 0.06 sin [27rt/10 V + 90°] 

obliquity of the state is given by Eq. (13). For the value of w chosen for the equations 
of motion determinations, the precession period would be 286,660 years and /' = 7.51°, 
whereas for the real value oi w — w^^ relative to the Laplace plane normal, the precession 
period would be 328,000 years and F = 8.6° (Yseboodt and Margot, 2005). Both sets yield 
wsin/' = 1.641 X 10-^°/ycar in the numerator of Eq. (13). With C/MR^ = 0.34, a central 
value of a plausible range (Harder and Schubert, 2001), e = 0.206, n = In j (87.969 d), 365.2563 
d/y, (J2, C22) = (6.0 X 10~^ 1.0 X 10"^) (Anderson et al. 1987), wR = 0.1408°/year and 
wS = 0.01437°/year (Recall R' and S have w in their denominators.), the two choices of w and 
corresponding /' yield ic = 1.6704' and 1.6696' respectively, a difference of only 0.05", which 
justifies the use of w in place of in the determinations of the Cassini state. Yseboodt and 
Margot (2005) find ic — 1.68'. If we use the published uncertainties in J2 and C22 of ±2 x 10~^ 
and ±0.5 x 10^^ respectively and assume that the uncertainties in the ephemerides used to 
determine /' and the uncertainty in C/MF^ are neghgible, 1.18 <ic< 2.51 arcmin. 

The unit vector corresponding to the Cassini state position ec is in the plane determined 
by the orbit normal and the e^^ on the opposite side of from e^. From the constraints 
that Be • (bo X e^) = and Gc • Go — cosic and the unit magnitude of Be, three equations in 
the three components of Be can be solved for the current position of the Cassini state. It is 
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Fig. 9. — Variation of e, I and Q, relative to the ecliptic of J2000 for JPL Ephemeris DE 408. 

less algebraically taxing to solve for the increments AX and Ay relative the projection of e,, 
on the ecliptic plane to yield the ecliptic latitude and longitude of the Cassini state position 
of Ac = 82.9694° and 0^ = -41.7585° respectively for the J2000 epoch with = 1.67'. The 
direction of the displacement is found from the intersection of the plane determined by 
and any of the w from Eq. (7) and the unit sphere. This position will of course have to be 
adjusted as the uncertainties in the determining parameters are reduced. 

The response of the angular separation 5 of spin axis and the Cassini state to short 
period fluctuations can be investigated by flxing e = 0.19 and dQ/dt = —27r/(287, 000 years) 
but forcing / to vary as / = Jq + sin (27rt/Pf). We fix the amplitude Aj, and determine 
the maximum separation of the spin from the Cassini state as a function of the period of 
the variation Pj. Panel d in Fig. 7 is characteristic of the behavior of the spin-Cassini state 
separation for periodic variations in e and /. The separation starting at 5o — fluctuates 
between zero and a maximum, which is about 1.2" in panel d in Fig. 7. Fig. 10 shows this 
maximum separation for periods ranging from 200 years to 5 x 10^ years for several values of 
Aj. Quite large separations of the spin from the Cassini state can result from such relatively 
short period oscillations in /, but most noticeably at periods near the period of the spin 
precession about the Cassini state, where a clear resonant response is evident. 

We saw in Section 4 how a very small spectral power in the Quinn data at a period near 
the resonance more than doubled the maximum spin-Cassini state separation over the 3 mil- 
lion year interval from the 1" maximum separation when spectral power at that period was 
suppressed by halving the number of data points. To check whether there is significant power 
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Fig. 10. — Maximum separation of the spin axis from the Cassini state for periodic variation 
of / with amphtude Aj with e = 0.19 and dfl/dt = —2tt/ (287, 000 years) as a function of the 
period of variation. 

in the real variations in e, /, Q, dl/dt and dVt/dt at the resonant period or at any other short 
period that could cause the spin to separate from the Cassini state by a large angle, we repeat 
the exercise in Section 4, but now with the variations in e, / and Q and their derivatives given 
by the ephemeris DE 408. Fig. 9 shows that the short period fluctuations in e, / and Q are 
small compared to significant almost linear variations. 

As was done with the Quinn data, we represent the variations in e, /, and with spline fits 
that capture the short period fluctuations in these elements as well as their derivatives derived 
therein for use in Eqs. (10). The Cassini state position is determined at arbitrary times from 
averaged values of the parameters, with (e), (/), and (Q) being averages determined by the 
sum of the values at the extremes of a 2000 year window centered on the epoch divided by 
2, and (dl/dt) and {dfl/dt) are the difference in the values at the extremes divided by 2000 
years. 

In Fig. 11, the equivalent of Fig. 5 the projections of the components of and ec onto the 
orbit plane are shown over the 20,000 year interval of the ephemeris, where the projection of 
Cc is truncated at the endpoints because of limits of the averaging process. Fig. 11 represents 
a short segment of one of the loops in Fig. 5 extended into the future to JIOOOO, except now 
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Fig. 11. — Trajectories of the projections of the unit spin vector and the unit vector in 
the Cassini state direction ec on the orbit X'Y' plane for variations in the orbital elements 
according to the 20000 year JPL Ephemeris DE 408. 

all of the high frequency terms are included in the variations of e, I and Q for the variation 
of the spin position. Again the spin is initially in the Cassini state {6o = 0), and we see that 
it remains close to that state as we have defined it in terms of the averages of the ephemeris 
data. In Fig. 12, the equivalent of Fig. 6b, we show that the spin, initially in the Cassini 
state, stays within 1" of the Cassini state over the 20,000 year interval of the ephemeris. There 
are no short period terms in the real variations of e and I that can lead to the large periodic 
separations shown in Fig. 10. In particular, there are no contributions to the variations that 
are near resonance with the spin precession, and the spin will remain close to the Cassini state 
for both long period and short period variations in the orbital elements and variations in the 
Laplace plane orientation that define the position of the state. 



6. Discussion 



We have shown that Mercury's spin axis will stay within approximately 1" of the time varying 
position of Cassini state 1 after dissipative processes have brought it to the state for both long 
and short period variations in the orbital parameters. Yseboodt and Margot (2005) also find 
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Fig. 12. — Separation 6 of the spin vector from Cassini state 1 for variations in e, / and 
from JPL ephemeris DE 408 covering 20,000 years centered on calendar year and sampled 
every 500 ephemeris days. The spin is initially at Cassini state 1 in this example, and the 
separation is deduced from a spline fit through the 14,624 data points of the ephemeris. The 
periodicity in S is that of the free precession of the spin about the Cassini state. 

that Mercury's spin axis remains close to the Cassini state as the position of the state is altered 
by slowly increasing the planetary masses from zero to their current values with time scales 
that are long compared to the spin precession period. The maximum separation of 1" that we 
found is about 1% of the 1.67' current obliquity of the Cassini state. Current estimates for 
the precision of the MESSENGER spacecraft determination of Mercury's obliquity are about 
10% of its value (Zuber and Smith, 1997), although this precision may possibly be exceeded 
by radar measurements. In any case the maximum deviation of Mercury's obliquity from the 
Cassini state obliquity induced by the slow and fast orbital variations and by the changing 
geometry of the solar system is sufficiently small that reasonably precise values of C/MR^ 
will be obtainable from the radar and spacecraft observations with correspondingly tight 
constraints on Cj^/C. This conclusion depends of course on there being no recent excitation 
of a free spin precession. With tight constraints on both C/ME? and C^/C, the internal 
structure of Mercury should be reasonably well constrained. 

For the small coupling that is likely between Mercury's hquid core and solid mantle, the core 
is not likely to follow the mantle on the spin precession time scale as we have assumed here. 
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If the mantle precesses independently of the core, the polar moment of inertia for the spin 
precession will be about half of the total moment of inertia. The spin precession period will 
be 500 years instead of the 1000 years assumed here. The mantle would relax to a somewhat 
different Cassini state which would gradually relax to the current state as dissipative processes 
cause the liquid core to catch up on a time scale short compared to the 300,000 year orbit 
precession time scale. It may be the case that there is a slight offset of the spin axis from 
the Cassini state because of the fluid coupling between core and mantle. If this offset is 
measurable, it will provide a constraint on the core-mantle coupling. Investigation of this 
conjecture will be included in a following paper, where core and mantle are considered as 
independent but coupled entities. 
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